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ABSTRACT 

We examine the effects of gas expulsion on initially sub-structured and out-of- 
equilibrium star clusters. We perform iV-body simulations of the evolution of star 
clusters in a static background potential before adjusting that potential to model gas 
expulsion. We investigate the impact of varying the rate at which the gas is removed, 
and the instant at which gas removal begins. 

Reducing the rate at which the gas is expelled results in an increase in cluster 
survival. Quantitatively, this dependency is approximately in agreement with previous 
studies, despite their use of smooth, and virialised initial stellar distributions. 

However, the instant at which gas expulsion occurs is found to have a strong effect 
on cluster response to gas removal. We find if gas expulsion occurs prior to one crossing 
time, cluster response is poorly described by any global parameters. Furthermore in 
real clusters the instant of gas expulsion is poorly constrained. Therefore our results 
emphasis the highly stochastic and variable response of star clusters to gas expulsion. 
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1 INTRODUCTION 

Many young stars are found in clusters. The fraction which 
are found in clusters depends on ones definition of a clus- 
ter (see Lada & Lada 2003; Bressert et al. 2010), and how 
exactly a cluster is defined (b ound vs. unbound, see e.g. 
iGieles fc Portegies Zwartl l201ll h However one looks at it 
though, star 'clusters' are clearly an important mode of star 
formation, if not the dominant mode. 

Most stars appear to form in groups of tens to thou- 
sands of members embedded in molecular clouds. The initial 
distribution of stars follows the complex clump y and fila- 
ment ar y structure of the u nderly in g gas (see e.g. Kirk et all 



ment ary structure o r tne underlying gas (sec e.g. rvirk et_ay 
20071: lOutermuth et ail 12009: fPeretto fe Fuller! 120091; 



Bressert et al.l 2010l : di Francesco et ail 2010l : Maury et al.l 
201ll ). Simulations of initially smooth clouds of molec- 
ular gas, that are seeded with su personic turbulence, 



fragment into filamen t ary s t ructu r es |Bonnell et al. 120031 : 
iBate fe BonneiH |2004 iBatd 120091 : iGirichidis et all I2OI1I ). 



These structures may collapse to form stars if they become 
self-gravitating. It has been argued that clusters form as 
entities within this complex distribution (e.g. Marks & 
Kroupa 2011), or that they can form afterwards by mergers 
of stars and stellar groups from this complex distribution 



I Allison et al.ll2009l : [Allison et al.l 120101 ') . Clearly the initial 
dynamical state of the young stars will play a crucial 
role in forming bound clusters or unbound asso ciations 
|Allison et al.ll2009l : IGieles fe Portegies Zwart]|201ll ) . 

In this paper we will take the view that clusters form 
within the hierarch ical merging scenario (|Allison et al 1 120091 : 
lAllison et aHl2010l ) and, together with gas not used in star 
formation, they form bound objects. Observations show 
very few young clusters associated with natal gas older 
than ~ 5 Myr jLada fc Ladal 120031 : lLadal I2010T I. It is as- 
sumed that feedback from massive stars removes residual 
gas in a gas expulsion phase. This gas expulsion will sig- 
nificantly alter the potential felt by the st ars and can re- 
sult in the destruction of the cluster (e.g. iTutukovl [l978l: 
Hills|ll980l:lGoodwin fc Bastianll200t3 : lBaumgardt fc Kroupal 
20071 : iGoodwinl l2009l h" This is often cited as the cause of 
'infant mortality': t he apparently high destruction rate of 
young clusters^] (e.g. lLada fc Ladall2003T ). 

Previous work on gas expulsion has tended to concen- 
trate on clusters in which the g as and stars are in virial and 
dynamical equi librium (but see I Verschueren fc Davidlll989l ; 
I Goodwinl [20091 ). The assumption of virial and dynamical 
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1 Although it should be noted that the importance of infant mor- 
tality depends on one's definition of a cluster. 
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equilibrium for young clusters can be used to infer the initial 
cluster pop ulation and properties from present-day popu- 
lations fe.g. IParmentier fc Gilmord 120071 : iBaumgardt et aD 
120081 ; IParmentier et alT 20081 ). However, the assumption of 
equilibrium initial conditions means that there is a one-to- 
one correlation between survivability and the initial condi- 
tions. 

This paper is part of a series in which we examine gas 
expulsion in the context of the hierarchical merging scenario 
in which the initial distributions of the stars and gas are not 
in dynamical equilibrium. Specifically they have an initial 
distribution which will cause them to collapse and form a 
bound cluster with in the gas potential (|Allison et al J 120091 ; 
lAllison et ai]|2010l : Smith et al. 2011a,b). 

In this paper we present a simple numerical experiment 
examining the evolution of non-equilibrium small- TV clusters 
within a smooth background potential which models the gas. 
Whilst this is clearly not realistic (at birth stars follow a gas 
distribution which is extremely complex, see above) we wish 
to examine the stochasticity which complex stellar distribu- 
tions introduce. As we shall see, statistically identical non- 
equilibrium clusters can evolve in very different ways even in 
our very simple numerical experimental and this stochastic- 
ity removes the simple correlation between initial conditions 
and survivability that is so often assumed. 

In this paper we examine the A r -body evolution of 
highly non-equilibrium star clusters in background poten- 
tials designed to model the residual gas in young star clus- 
ters. Here we adjust the background potential to model the 
effects of gas expulsion on different timescales, starting at 
different times. We describe our initial conditions in Sec- 
tion 2, our results in Section 3, and then discuss the potential 
consequences in Section 4 before drawing our conclusions in 
Section 5. 



2 INITIAL CONDITIONS 

In this paper we present a very simple numerical experi- 
ment: a clumpy distribution of equal-mass stars moving in 
a smooth (Plummer) background potential which is then 
removed to simulate gas expulsion. This experiment is not 
meant to accurately reflect reality, rather it is meant to con- 
centrate on a new important variable: complex initial stellar 
distributions which can introduce (very significant) stochas- 
ticity to the problem of gas expulsion. As we will discuss 
in Section 4 we feel it emphasises some crucial physical pa- 
rameters that would be important in any more sophisticated 
simulations and, most importantly, in reality. 

We perform our Af-body simulations using the Nbody6 
code (|Aarsethll2003l ). In our simulations there are two sep- 
arate mass distributions which we wish to model: the stars 
and the background gas potential. 



lar mass of 5OOM0 . We choose equal- mass particles in order 
to av oid complex two-bo dy in teractions and mass segrega- 
tion. lAllison et all (120091 ) and lAllison et all (|2010l ) showed 
that both of these effects can be extremely important in the 
violent collapse of cool, clumpy regions. However we wish to 
avoid complicating our simulations with these effects as we 
are interested in the effects of gas expulsion. 

We distribute the stars within a radius of 1.5 pc with 
a fractal distribution. We choose the fractal dimension D to 
be D — 1.6. The fractal is constructed using the box fractal 
method described in detail in Goodwin & Whitworth 2004 
(see also lAllison et al.ll20l61 ). A fractal with D = 1.6 corre- 
sponds to a highly clumpy initial distribution (see upper-left 
panel of Figure [2]). It should be noted that clumpy fractal 
clusters can vary considerably in appearance depending on 
the random realisation used a nd their subsequent evolution 
can be highly stochastic (see lAllison et al. I l2010l ; Smith et 
al. 2011a). When investigating their response to varying gas 
expulsion time-scales, we therefore conduct a minimum of 5 
realisations of any cluster. 

The gas within the star-forming region is modelled as a 
background Plummer potential with a mass M g = 3450 Mq 
and scale-radius r g = 1.0 pc. The crossing-time of the region 
is t cr ~ 1.3 Myr. 

We emphasise that the gas potential does not follow 
the initial stellar distribution, nor does it react to changes 
in the stellar distribution as it evolves (it is not live). These 
are obviously extreme simplifications but, as we will discuss 
later, we feel that we capture the essence of the basic physics 
using such a simple model. 

We set the initial velocity dispersion of the stars rela- 
tive to the total potential (gas & stars) with initial stellar 
virial ratios, Qj, of 0.2 (sub-virial), or 0.5 (virialised). Note 
that a clumpy stellar distribution with Qi = 0.5 is not ini- 
tially in dynamical equilibrium. Thus the cluster will evolve 
significantly due to violent relaxation. 



2.2 Star formation efficiency 

The choice of the initial mass (M g ) and Plummer radius (r g ) 
of the (Plummer) gas potential sets the true star formation 
efficiency (SFE), i.e. initial fraction of gas that has been 
turned into stars. For a gas mass of M g = 3450 Mq, and 
scale radius r s — 1.0 pc, within 1.5 pc the enclosed gas mass 
is 2000 Mq. The stars are distributed out to 1.5 pc, with a 
total stellar mass of M s = 500 Mq . Thus within a radius of 
1.5 pc the true SFE is 20 per cent. 

We note that this is below the critical SFE for the sur- 
vival of a bound core from an initially virialised star-gas dis- 
tribution with instantaneous gas expulsion (~ 30 per cent, 
iGoodwin fc Bastianll2006l ; IBaumgardt fc Kroupall2007l ). and 
at the lower bound f or the survival of a bound co re for adi- 
abatic gas expulsion (|Baumgardt fc Kroupall2007l ). 



2.1 Initial distributions 

In all cases we model the stellar distribution as N = 1000 
particles with equal masses of 0.5Mq resulting in a total stel- 

2 And it is difficult to see how the addition of more realistic 
physics can make the problem less stochastic rather than more. 



2.3 Gas expulsion 

Although young star clusters form embedded within the 
molecular gas from which they formed, few sta r clusters over 
~ 5 Myr old remain associated with their gas (|Lada fc Ladal 
120031 ; iProszkow fc Adamsll2009l ). This is likely as a result of 
a number of mechanisms including radiative feedback from 
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Figure 1. Plot of evolution of gas mass with time using Equa- 
tion [l] (black curves), assuming i cr =l-3 Myr. td ur is the length 
of time required for the gas mass to be reduced to 99% of its 
initial value. When td ur is short, gas loss occurs over a short 
duration. Our slowest gas mass-loss is when t^ UT =50 Myr. For 
comparison, we includ e curve s of the evolution of gas mass from 
IBaumg ardt & Kroupa 1 200 3) (red curves). T m is the parameter 
from lBatuTigard*t"*irKr"oupal l|2007l ) controlling the rate of gas loss. 
We see that td ur =50 Myr results in similar gas mass-loss to when 
r m =l-2 crossing-times. 



massive stars, stellar winds from young stars, and eventu- 
ally the first supernova(e). The time at which gas removal 
begins to occur, and the duration of the gas removal process 
is uncertain, and dependent on the particular gas removal 
mechanism in operation. 

One of the main aims of our simulations is to examine 
the decoupling of the gas and stars before gas expulsion. If 
the gas and/or stars are not in equilibrium either as a whole, 
or with each other, then the relative distributions of the 
gas and stars will change. As shown in our previous papers 
(Smith et al. 2011a,b) if the stars can collapse relative to 
the gas distributio n then the effec t of gas expulsion is less 
important (see also lGoodwin|[2009l and references therein). 

To simulate gas expulsion we parameterise the time evo- 
lution of the mass of the gas potential as: 



M g (t) 



M g (0) 



1 + M(t - t e 



(1) 



where t cxp is time at which gas expulsion begins, and M 
controls the rate at which gas is removed. 

We choose M such that the gas mass M s (t) is reduced 
to 1 % of its initial value M g (0) over a specified duration 
time tdur- For example, if t axp = 1 Myr and tdur = 5 Myr, 
then gas expulsion begins 1 Myr after the beginning of the 
simulation, and 99% the gas will be gone 5 Myr later. 

We test the response of a cluster to four values when 
gas expulsion begins (i exp : 0.1, 1, 5 & 9 Myr), and four val- 
ues for how fast gas expulsion occurs (fd ur : 0.01, 1.0, 5.0 
& 50.0 Myr). These timescales of gas expulsion cover the 
reasonable range of possible gas expulsion timescales from 
instantaneous to extremely slow (adiabatic) compared to 



the crossing time. It is unclear what the typical gas expul- 
sion timescales are, an d they may depend on cluster mass 
(|Kroupa fc Boily||2002h . For example, the relatively late su- 
pernovae of a low-mass O-star might result in a delayed start 
to gas expulsion. 

In Figure [TJ we illustrate the reduction in gas mass for 
gas expulsion that begins effectively immediately (t cxp : 0.1 
Myr), and occurs over a duration td u r : 0.01, 1.0, 5.0 fc 50.0 
Myr (black curves). We note that IBaumgardt fc Kroupa! 
(|2007l ) use an exponential form for the evolution of the gas 
mass tha t differs from ours in de t ail, b ut not in essence. 

The IBaumgardt fc Kroupa! ([2007) form is shown on 
Figure [1] as red curves (labelled BK07). We see that our 
slowest gas mass reduction (id ur =50 Myr) can be con- 
sidered approximately equiv alent to the gas mass loss in 
IBaumgardt fc Kroupal (|2007l ) when r m =l-2 crossing-times. 
We emphasise that while a 50 Myr duration (i.e. until 99 
per cent gas loss) for gas expulsion may sound unrealisti- 
cally lengthy, this is merely a result of the analytical form 
that we have used to model gas expulsion. 



2.4 Cluster properties at time of gas expulsion 

For a star-gas cluster which is initially smooth and in virial 
equilibrium (i.e. in dynamical equilibrium) and is no longer 
forming stars, its response to gas expulsion can be char- 
acterised by two parameters. Firstly, the true SFE, the 
relative gas mass to stellar mass (which does not change 
as the cluster is in equilibrium). Secondly, the timescale 
of gas expulsion from instantaneous (in less than a cross- 
ing time), to adiabatic (in s everal crossing times). See e.g. 
IBaumgardt fc Kroupal (I2007T) 

However jGoodwir] (|2009l ) points out that the important 
variables are actually the virial state of the stars at the onset 
of gas expulsion and the timescale of gas expulsion. In an 
initially virialised cluster the virial state of the stars can be 
found directly from the true SFE. However if the cluster 
is not in equilibrium initially, then the virial state of the 
stars corresponds to an effective SF E which can be very 
different from the true SFE (se e also IVerschueren fc Davidl 
ll989l ; lGoodwin fc Bastianll2006l ). 

In Smith et al. (2011a), we parameterise the effective 
SFE using two more physically meaningful quantities. The 
Local Stellar Fraction (LSF) is a measure of the mass of gas 
in the region where the stars are found (a form of evolving 
effective SFE). We also define the virial ratio of the stars at 
the moment of gas expulsion Q ge . 

The LSF is defined as 



LgF _ M *(r < r Hs) ) 

AL(r < r h(s) ) + M g (r < r h ( H )) 



(2) 



and 



where r h ( s ) is the half-mass radius of the stars, and M* 
M g are the mass of stars and gas, respectively. 

The key aspects of the LSF and Q gc are that they evolve 
with the dynamical evolution of the cluster. In Smith et al. 
(2011a) we found the LSF and Q ge (measured when gas 
expulsion begins) to be effective predictors of the response 
of an initially out-of-equilibrium cluster to instantaneous gas 
expulsion after a few crossing times (in their case 2.5 initial 
crossing times). 

We measure the survival of clusters to gas expulsion by 
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Figure 2. The evolution of an initially fractal (D = 1.6) stellar distribution in an r g = 1 pc, M g = 3450Mq gas potential. The initial 
virial ratio of the stellar distribution within the gas background potential is Qi=0.2. The time of each snapshot is indicated in the top- 
right hand corner of each panel. Initial clumpy substructure is erased with time. Instantaneous gas expulsion occurs at 5 Myr, resulting 
in unbound stars. By 15 Myr the final cluster is settled. 
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Figure 3. Final bound fraction of stars /bound versus Local Stel- 
lar Fraction (LSF) for clusters where t e xp=4-7 t CT . Symbol colour 
denotes the duration over which gas expulsion occurs (see legend). 
A trend for increasing /bound with increasing LSF can be seen. 
If gas expulsion occurs more slowly, this trend shifts vertically 
upwards on the diagram. By eye, our slowest gas removal rate 
appears to result in an approximate increase in /bound of ~ 0.3. 



measuring the fraction of stars that are bound to the clus- 
ter at t = 15 Myr. Where sub-structure remains following 
gas expulsion, we choose the most massive sub-clump to be 
the main cluster. We measure the number of stars that are 
bound to the cluster within the inertial frame of the cluster. 
This is because some clusters have a net velocity though 
simulation space following gas expulsion. To calculate the 
average velocity, we use only stars within a specific radius 
of the cluster centre. We vary the radius to ensure our fi- 
nal bound fraction is not sensitive to the radius we choose. 
When measuring the bound fraction we only include stars 
that are bound to the main cluster, and exclude stars that 
are bound to small subclumps in any remaining substruc- 
ture. 



2.5 Summary 

We set-up clumpy 5OOM0 star clusters with N = 1000 equal- 
mass particles within a static background gas potential. The 
initial dynamical state of the stars varies from cool and sub- 
virial (Qj=0.2) to virialised (Qi=0.5). We then allow the 
stars to dynamically evolve within the gas potential. Gas ex- 
pulsion begins at a time t — t oxp . We vary t exp from 0.1 Myr 
(essentially immediate gas expulsion) to 9 Myr (~ 7 crossing 
times). At the start of gas expulsion, we measure the Local 
Stellar Fraction of the cluster (LSF), and stellar virial ratio 
(Q gc ). We also vary the rate at which gas expulsion occurs 
by setting the parameter t^ur (the amount of time required 
for 99% of the gas to be expelled). We consider 4 values: 
Wr = 0.01, 1.0, 5.0, & 50 Myr. We measure the final bound 
fraction of each star cluster at t = 15 Myr to determine how 
well the cluster survives. 




2.5 3 3.5 4 4.5 

Crossing-times (t/t ) 

Figure 4. Evolution of the stellar virial ratio of 5 different fractal 
clusters. All start with a total virial ratio of Q; = 0.5. 



3 RESULTS 

3.1 A representative example 

We will first present a representative example, before fo- 
cussing on the impact of the gas expulsion times scales. In 
Figure [5] we present an example of the evolution of an ini- 
tially fractal D = 1.6 stellar distribution with an initial virial 
ratio of Q\ — 0.2 in the gas potential. The true SFE of this 
system is e = 0.2. In this particular example gas expulsion 
begins at 5 Myr, and then the gas is removed effectively in- 
stantaneously. The crossing-time of the region is ~ 1.3 Myr. 

Figure[2]shows that the stellar component collapses and 
relaxes during the first two crossing times. The initial fractal 
distribution (first panel) is partly erased by one crossing 
time (second panel) and by two crossing times (third panel) 
the bulk of the clu ster has collapsed into a dense core (cf. 
lAUison et alll2009l ; Smith et al. 2011a,b). At the onset of gas 
expulsion at 4.5 Myr (fourth panel) the stellar distribution 
is relatively relaxed with a halo of stars, mostly ejected by 
two-body interactions within the dense main stellar cluster. 
The effect of gas expulsion is as expected in that a significant 
fraction of stars are lost. At 5.9 Myr (fifth panel) many of 
these escaping stars are still associated with the cluster, but 
by 15 Myr (last panel) only a virialised bound core of stars 
remains. 

The cool, clumpy initial distribution of the stars in this 
simulation cause them to collapse towards the centre of the 
gas potential. This causes the stars to virialise (i.e. Q &e to 
approach 0.5), but it also causes the LSF to increase signifi- 
cantly as the stars are concentrated in the centre of the gas 
potential. For example, in Figure [2] the initial distribution 
of stars has an LSF=0.24 (close to the true SFE), but the 
final LSF at the time of gas expulsion has increased to 0.41. 



3.2 



Effects of varying the rate of gas expulsion 



In all simulations, we evolve our stellar initial conditions in 
the gas potential until gas expulsion begins. At the start of 
gas expulsion we measure the LSF (see Equation [5} and the 
virial ratio Q gc . 

We simulate the evolution of 10 different realisations. 
For each ensemble we model gas expulsion occurring at t cxp = 
0.1, 1.0, 5.0 & 9.0 Myr. For each chosen i C x P , we model four 
durations over which the gas is expelled: tdur= 0.01, 1.0, 
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5.0, & 50 Myr. We use Qi=0.2 and 0.5. In total, we conduct 
90 simulations. All simulations are conducted for 15 Myr in 
total (~ 11.5 initial crossing-times), when the final bound 
fraction /bound of the cluster is measured. 

The results of the t exp =5. 0-9.0 Myr (£ C x P =4-7 crossing- 
times) simulations are shown on Figure [3] (we turn later to 
the effects of different t cxp ). We plot the final bound frac- 
tion /bound of the cluster against the LSF at the start of 
gas expulsion. The colour of the filled circle symbols de- 
notes the duration over which their gas was removed. For 
clarity in Figure [3] , we exclude the tdui = 1 Myr points, 
this gas removal timescale is effectively instantaneous and 
the results are indistinguishable from the tdur = 0.01 Myr 
points. Black symbols have almost instantaneous gas expul- 
sion (/dur=0.01 Myr). Red symbols lose their gas marginally 
more slowly (idur=5 Myr), and green symbols lose their gas 
the most slowly of all (£dur=50 Myr). (For reference see Fig. 

m) 

A trend can be seen in the relation between /bound and 
the LSF in Figure [3] Firstly, as expected, lower LSFs tend 
to result in lower bound fractions. Also, as expected, slower 
gas removal (tendin g to adiabatic) tends to resu lt in higher 
bound fractions (cf. IBaumgardt fcK roupa 2 0071 ). 

However, for low-LSFs there is a wide spread in the final 
bound fractions. For LSFs of ~ 0.2 (similar to the initial, 
true, SFE), the bound fractions vary between zero and over 
half. The largest bound fractions (all > 0.4) occur for adia- 
batic gas loss which is far less destructive (see above), but 
for an intermediate gas loss timescale of 5 Myr the spread 
is between zero and 0.3. We will return to this in the next 
section. 

We do find broad ag reement with the results of 
IBaumgardt fe Kroupal (|2007l ) once the differences in our pa- 
rameterisation of the gas expulsion timescale and our use of 
the LSF rather than true SFE are taken into account. Our 
LSF values are typic ally between 0.2 and . 4. W e there- 
fore compare with the IBaumgardt fc Kroupal (|2007l ) results 
for true SF Es in the same range. The u pper-left panel of 
Figure 2 of IBaumgardt fc Kroupal |2007l ) shows that their 
bound fraction also increases by 20-40 per cent between their 
effectively instantaneous and their adiabatic gas expulsion 
timescales. Inspection of Figure [3] shows a similar relation- 
ship, albeit with much more scatter. 



3.3 The evolution of the dynamical state of the 
stars 

The scatter we find in the LSF- /bound relationship for a given 
gas expulsion onset and duration comes from two sources. 
It is due to differences in the dynamical state of the stars 
at the onset of gas expulsion which, in turn, is due to the 
stochastic evolution of the initially non-equilibrium stellar 
distribution. 

Figure U shows the evolution of the stellar virial ratios 
for 5 different fractals that initially have a total virial ratio of 
Qi = 0.5. The stellar virial ratio is the virial ratio of the stars 
alone. We note that the potential of any individual star has a 
component from surrounding stars and another component 
from the surrounding gas potential. Both components are 
included in the virial ratio calculation. 

Although the systems are initially in virial equilibrium, 
the fractal initial conditions mean that they are not in dy- 
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Crossing times (t/O 

Figure 5. Evolution of peak values of |(Q(t) — 0.5) | for 10 dif- 
ferent initially virialised (Q;=0.5) fractals and 10 different cool 
(Qi=0.2) fractals. The average (and standard deviation) of the 
peak values are binned into half crossing-time time-bins. Peak 
values arc initially larger, then reduce in size as the clusters at- 
tempt to virialise. Initially Qi=0.2 fractals have larger deviations 
from virialised. After ~ 1.5 crossing-times, Q-cvolution of Q;=0.2 
& 0.5 are similar. After 4-5 crossing-times, fractals have settled 
to close to virialised. 



namical equilibrium. The stellar distributions violently relax 
and collapse which initially increases the stellar virial ratio 
before finding a roughly smooth, dense new configuration 
in the centre of the gas potential (with a correspondingly 
higher LSF). The stars then 'bounce', oscillating around 
virial equilibrium for several crossing times whilst they re- 
lax. 

Smith et al. (2011a) found that the stellar virial ratio 
at the onset of gas expulsion is im portant to the c luster's 
response to gas expulsion (see also iGoodwinl 120091 ) . If the 
stellar virial ratio is slightly super- virial (Q go > 0.5), then 
the cluster will loose more stars (lower /bound), whilst if 
they are slightly sub- virial (Q ge < 0.5), then the cluster 
will retain more stars (higher /bound) - However, Smith et al. 
(2011a) only examined instantaneous gas expulsion after 2.5 
crossing times. 

To investigate the effects of Q go with different starting 
times and durations of gas expulsion we record the peak 
values of the deviation of the stellar distribution from virial 
equilibrium (|(Q — 0.5) | ) and the time at which these maxi- 
mum deviations occur. We conduct this test both for 10 dif- 
ferent clusters that are initially globally virialised (Qi=0.5), 
and 10 different clusters that are initially globally cool 
(Qi=0.2). We calculate the mean and standard deviation 
of the peak values of \(Q — 0.5) in bins of half an initial 
crossing-time. The results are presented in Figure [5] 

Figure [5] quantifies what is visible in Figure [4] Devia- 
tions from virialised are greater at early times as the stellar 
component relaxes. They are also greater in initially globally 
cool clusters (Qi=0.2) as they have to relax more to reach 
(rough) equilibrium. 

We can predict that the stellar virial ratio will be an 
important factor in that it will vary more, and be more 
stochastic, at early times. Firstly, the further clusters start 
away from equilibrium (dynamical or virial) , the greater the 
initial variations can be. Secondly, the shorter the amount 
of time clusters have had to relax, the greater the variations 
can be. 
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Figure 6. Final bound fraction of stars /bound versus Local Stellar Fraction (LSF) for clusters with early gas removal (t C xp < 1 ten 
upper panel), intermediate gas removal time (t CX p=1.5-2.0 tcr, middle panel), and late gas removal (t e xp=4-7 t C r, lower panel). In these 
models, when gas expulsion begins, it is effectively instantaneous (idur=0-01 Myr). Clusters with Qi=0.5 result in lower LSF values, 
whereas Q\=Q.2 clusters result in higher LSF values. The stellar virial ratio at the onset of gas expulsion Q gc is either low (blue circles, 
Q ge < 0.46), roughly virialised (open circles, 0.46 < Q ge < 0.54), or high (red squares, Q gc > 0.54). Light crosses show the results from 
Smith ct al. (2011a) for instantaneous gas expulsion at i exp = 2.5 t cr - for reference. 
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Figure 7. Final bound fraction of stars /bound versus Local Stel- 
lar Fraction (LSF) for clusters with instantaneous gas expulsion 
at t eX p = 2.5 t cr - Results from Smith et al. (2011a) where the 
star formation efficiency is 20%, and using an initial distribu- 
tion of stars that is 'Plummer-like' (black circles), and fractal 
(red triangles). New results are shown for clusters with matching 
properties to the fractals (red triangles), but a much higher star 
formation efficiency of 60% (green squares). The scatter in the 
bound fraction is not highly sensitive to the form of the initial 
clumpy stellar morphology we choose, nor to the star formation 
efficiency. 



3.4 The survival of clusters 

From the above results and discussions we can see that the 
response of a cluster to gas expulsion, and the possible sur- 
vival of a bound core as a remnant cluster, depend on several 
factors. 

The relative importance of the gas potential to the stars 
(as measured by the LSF) is clearly very important. This, 
however, depends on the way in which the stellar distri- 
bution has relaxed relative to the gas and depends on the 
initial distribution (in this case highly clumpy) and initial 
velocity dispersion (measured by the total initial virial ra- 
tio Qi) of the stars. It also depends on the exact dynamical 
state of the stars at the moment when gas expulsion occurs 
(Q gc ). And both the LSF and Q gc depends to some extent 
on the stochastic details of the relaxa tion of each particu- 
lar realisation (see lAllison et af1l2010| ). The reaction of the 
cluster to gas expulsion also depends on the timescale of gas 
expulsion, with long (adiabatic) timescales being much less 
destructive. 

In addition, as we will explore in this section, both the 
LSF and especially Q gc evolve rapidly, especially in the first 
1-2 crossing-times and so the effect of gas expulsion will 
depend strongly on the time of the onset of gas expulsion. 

We have already investigated the effect of varying the 
rate of gas expulsion in Section 13.21 We now focus on the 
effect of varying the instant at which gas expulsion occurs, 
and fix the rate of gas expulsion to be effectively instanta- 
neous. 

Figure [6] shows the LSF- /bound relationships for early 
gas expulsion (top panel, t exp < 1 t cr ), intermediate gas ex- 
pulsion (middle panel, t cxp = 1.5 — 2.0 t CI ), and late gas ex- 
pulsion (bottom panel, t cxp — 4 — 7 t CI ). For all data points 
shown, when gas expulsion begins, the gas is removed in- 
stantaneously. Plotted on each panel in light crosses are the 
results from Smith et al. (2011a) instantaneous gas expul- 



for Qi=0.5 and Qi=0.2 star clusters in Figure [6] 

The points marked with blue filled circles have a low 
stellar virial ratio at the moment of the onset of gas ex- 
pulsion (Qgc < 0.46), open circles are roughly virialised 
(0.46 < Qgc < 0.54), and red squares are super- virial 
(Qgc > 0.54). 

The top panel shows the results for the early onset of gas 
expulsion. Early gas expulsion results in gas expulsion oc- 
curring whilst the stars are still relaxing. They tend to have 
a low virial ratio (Q ge < 0.46) as they are still collapsing 
into a new equilibrium, and their LSF tends to be low (close 
to the true SFE) as they are still fairly widely distributed in 
the gas potential. The huge scatter in this top panel is due 
to the highly stochastic nature of the early dynamics of the 
clusters. Many of these clusters will still contain subclumps 
(see Figure [2]) and the bulk dynamics of these subclumps is 
important. At these stages Q gc is not a good indicator of 
survival as much of the stellar kinetic energy can be con- 
tained in bound clumps. If these clumps escape due to gas 
expulsion /bound will be low (or even zero), however retain- 
ing one or two subclumps can result in a very high /bound- 
Clumps also raise the stellar-to-gas fraction very locally at 
their centres. The LSF is insensitive to this as it measures 
the stellar-to-gas fraction on scales that are typically much 
larger than individual clumps (i.e. the half-mass radius of 
the stars). 

For intermediate gas expulsion timescales (middle 
panel) the LSF becomes more important as the clusters tend 
to have had the time to relax into a smooth central cluster. 
More clusters are roughly virialised and there is a trend of 
/bound decreasing with lower LSF as would be expected. As 
clusters are oscillating around virial equilibrium within the 
gas potential, the spread in /bound for a given LSF is wide. 
At a fixed LSF, sub-virial clusters tend to higher /bound 
than super-virial clusters. Therefore at intermediate gas ex- 
pulsion timescales, a combination of Q gc and LSF determine 
the response to gas expulsion. 

For long gas expulsion timescales (bottom panel) the 
clusters have had a chance to virialise within the gas po- 
tential and so the relationship between the LSF and /bound 
becomes much tighter with the LSF being the key factor in 
determining the response to gas expulsion. 

3.5 Dependency on initial stellar morphology and 
star formation efficiency 

First, we wish to test how sensitive our results are to the 
mathematical form of the initial clumpy distribution of stars 
we have chosen. To do this, we use models from Smith et 
al. (2011a) where we considered two types of clumpy dis- 
tributions (or stellar morphologies) for the initial positions 
of our stars; fractal (like those modelled in this paper), and 
'Plummer-like'. In the 'Plummer-like' distribution there are 
16 clumps of stars in total, and each clump is modelled as 
an individual Plummer sphere of mass ~ 30 M@, and Plum- 
mer scale-length ~ 0.05 pc. The clumps themselves are dis- 
tributed according to a Plummer distribution with a Plum- 
mer scale-length of 1-1.5 pc). Further details can be found 
in Smith et al. (2011a). 

In Figure [7] we show the final bound fraction of stars 
/bound versus Local Stellar Fraction (LSF) for clusters from 
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Smith et al. (2011a), but split into 'Plummer-like' (black cir- 
cles), and fractal (red triangles). We recall that both groups 
have equal star formation efficiency of 20%, and both un- 
dergo instantaneous gas expulsion after 2.5 crossing-times 
- they differ only in their initial stellar morphology. We 
find that the trend, and its scatter in bound fraction, is not 
highly sensitive to which form of clumpy stellar morphology 
we choose. Two apparently separated groups of Plummer- 
like models are seen at high and low LSF. However this is 
mainly just a result of low numbers of Plummer-like models 
with Qi =0.2-0.4, whose presence would fill the gap between 
the two groups. 

Secondly, we wish to test if our results differ when we 
consider clusters whose mass is dominated by stars, as we 
have so far only considered clusters that are heavily domi- 
nated by gas. Therefore we model 20 clusters with a much 
higher star formation efficiency of 60% (green squares). 
These differ from the fractals on Figure [7] in their star for- 
mation efficiency only. As expected, these clusters fall in 
the high LSF range. More importantly, we find that there 
is still a similar amount of scatter at higher star formation 
efficiency - despite the fact that the stars now dominate the 
mass of the cluster. In fact, the virial state at the time of gas 
expulsion Q ge continues to play the same role, broadening 
the trend, as when the star formation efficiency was 20%. 

In summary, the trend and its scatter is not very sensi- 
tive to the form of the initial clumpy stellar morphology we 
choose, nor to the star formation efficiency (i.e. the strength 
of the rigid background potential). 



4 DISCUSSION 

We have examined the response of young clusters to gas 
expulsion as measured by their ability to retain a bound core 
of stars and avoid complete destruction: their survivability. 

Clusters older than ~ 5 Myr are not seen to be as- 
sociated with the gas from which they formed. There are 
two possibilities for the fate of this gas. Firstly, that all of 
the gas was used in star formation (100 per cent efficiency). 
However this is unlikely as shortly after a massive (O- or 
B-star) forms it will heat the gas and remove it from the 
site of star formation (through winds, ionising radiation, or 
supernovae). The removal of a (very significant) fraction of 
the (gas) mass in a region will have a (very significant) effect 
on the potential felt by the stars that remain. If the stars 
are unbound or very loosely bound at formation they will 
be dispersed, and gas expulsion merely speeds this process 
somewhat. However, if the stars are bound then gas expul- 
sion can completely destroy the cluster, or leave a bound 
core which will be seen as a naked cluster. It must happen 
that bound cores are left after gas expulsion as many old, 
bound clusters are seen (most strikingly, globular clusters). 

A key question to address is what can be said about the 
initial conditions of old, bound clusters? What must be true 
is that they have survived gas expulsion and therefore they 
must have been bound at the onset of gas expulsion. 

Most previous work on gas expulsion has concentrated 
on clusters that are in virial and dynamical equilibrium at 
the onset of gas expulsion ( see e.g. ■Tutukov 1978; Hill j |l980l ; 
iGoodwin fc Bastianl [20061 ; IBaumgardt fc Kroupall2007f ). In 



a cluster is the true SFE, secondary factors are the timescale 
of gas expulsion (adiabatic gas expulsion is less destruc- 
tive), and the tidal radius (expansion over the tidal ra- 
dius can rapidly destroy a c luster), as shown in detail by 
IBaumgardt fc Kroupal (120071 ). 

However, observations of star formation lead us 
to believe that stars do not form in equilibrium. 
They certainly seem to form in a compl e x distribu- 
tion t hat follows the gas (e.g. iLarsonl 1 19951: iTesti et al 
2000l ; ICartwright fc Whitworthl |2004| ; IGutermuth et al 



20051 ; lAllen et al.l 120071; iKirk et all 120071; ISchmeia et al 
20081 ; IGutermuth et all 120091; IPeretto fc Fuller! 120091: 



Bressert et alj 2010l ; dT*Francesco et al.l 2010l ; iMaurv et all 
201ll ). they seem to f orm w i th a l ower velocity dis- 



persion than the 



Adams et all 120061 



the gas (lAndrel |2002|; [Walsh et all 
20061; IPeretto et all 120061; IKirk et al] 



2004; 



2007 



Proszkow fc Adamsll2009h. and their ear l y evolution appears 
to be dramatic (e.g . iBastian et"al1l2008l ; lAllison et al]l2009l ; 



such a situation the primary factor controlling the survival of 



lAllison et ai]|201Cft . 

In this paper we have shown that if the initial stellar 
distribution is out-of-equilibrium then the effect of gas ex- 
pulsion can be very stochastic (see also Smith et al. 2011a,b). 
In particular, there are several key factors: 

1. The effective SFE (measured by us as the LSF): the rela- 
tive stellar and gas masses at the onset of gas e xpulsion (see 
also lVerschueren fc Davidlll989l ; lGoodwinll2009h . The higher 
the effective SFE, the greater the survivability. 

2. The stellar virial ratio (Q ge ): the virial state of the stars 
alone at the onset of gas expulsion. The lower the stellar 
virial ratio, the greater the survivability. 

3. The time of the onset of gas expulsion: this effects to what 
extent the stars and gas can relax and/or decouple and so 
changes the effective SFE and stellar virial ratio. From our 
simulations, the earlier the onset of gas expulsion the greater 
the stochasticity in survivability. 

4. The timescale of gas expulsion: instantaneous gas 
expulsion is more destruc tive than adiabatic (see also 
IBaumgardt fc KroupalbOOTT l. 

5. Tidal radius: whilst we have not simulated this, there is 
no reaso n to think it should not have the same effect as that 
found by IBaumgardt fc Kroupall2007l . 

The first four of these factors are all stochastic and de- 
pend on the details of the initial conditions. 

We reiterate here that our initial conditions are not real- 
istic. We evolve equal-mass stars in a static background po- 
tential which we then remove in a simplistic way. However, 
we find it impossible to see how more realism could possi- 
bly produce less stochasticity. Indeed, more realistic initial 
conditions should introduce vast quantities of variables that 
are potentially important. 

For example, a realistic clumpy, live gas potential could 
react in very different ways. Depending on the details of 
the gas it could couple well to the stars (i.e. lower the 
LSF by moving with the stars), or couple badly (by being 
accre ted or removed from subclumps, e.g. iKruiissen et al.l 
l2012h . Depending on the structure of the gas, feedback could 
be less efficient, pos sibly escaping through bubbles (e.g. 
iDale fc Bonnellll201ll ). changing the onset of significant gas 
expulsion and/or its timescale. In addition, realistic small- N 
subclumps would contain a mass spectrum of stars, include 
binaries, and have a wealth of possible iV-body interactions 
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that are extremely stochastic on their own (|AlIison et al.l 
l2010h . 

In summary, there is no reason to think that if one 
cluster survives and another does not that they had very 
different initial conditions. Two sets of statistically identical 
initial conditions can produce a cluster that retains most 
of its stars or one that is completely destroyed. The final 
outcome depends on a wealth of different variables (only 
a handful of which we have rather simplistically modelled 
here), as well as the details of the initial conditions. 



5 SUMMARY & CONCLUSIONS 

We perform JV-body simulations of sub-structured, non- 
equilibrium 25OOM0 clusters with 500Mq of stars in N = 
1000 equal-mass particles. A static background potential 
provides a gas mass of 2OOOM0 surrounding the stars. The 
initial true star formation efficiency is fixed at 20 per cent. 
Gas expulsion is modelled as a reduction in the mass of the 
background potential over time. We allow the instant that 
gas expulsion begins, and the duration over which gas ex- 
pulsion occurs, to vary. 

We use two parameters to characterise the state of the 
cluster at the moment of gas expulsion. The Local Stellar 
Fraction (LSF) is a measure of the relative contribution of 
gas and stars to the potential felt by the stars. The stellar 
virial ratio is a measure of how relaxed the stellar component 
is. 

Our key results may be summarised as follows. 

(i) In initially clumpy, and non-equilibrium clusters, a 
slower rate of gas loss increases the survival of the cluster 
(in keeping with previous results). 

(ii) If gas is expelled within an initial crossing time the 
effect of gas expulsion is completely stochastic and there are 
no good global indicators of its effects. 

(iii) If gas is expelled after a few crossing times then a 
combination of the LSF (or rather the efficiency of the cou- 
pling of the gas to the stars), and stellar virial ratio are 
required to parameterise cluster survival, although it is still 
rather stochastic. 

(iv) If gas is expelled after several crossing times the stel- 
lar component is generally relaxed and the LSF alone is a 
good indicator of the effects of gas expulsion. However, the 
LSF can substantially differ from the true star formation 
efficiency due to the decoupling of the gas and stars. 

Our simulations have been very simple numerical ex- 
periments. However we feel they do capture some important 
physics. 

Firstly, the coupling of the stars and gas is crucial to 
the survival of a cluster. If the stars are able to relax to a 
more concentrated distribution than the gas (or conversly, 
if the gas is able to expand relative to the stars), then the 
dramatic effects of gas expulsion are reduced. 

Stochasticity is important: the exact level of decoupling 
from the gas, and the details of the stellar phase space dis- 
tribution are also crucial. Gas expulsion whilst the stars are 
still 'clumpy' can completely destroy one cluster, whilst leav- 
ing another - statistically identical - cluster with most of its 
stars in a bound cluster. This depends on the bulk virial ratio 
of the stars, and also the details of the clump velocities. It is 



only after several initial crossing times - after a cluster has 
relaxed into a rough equilibrium - that stochasticity reduces 
in importance. However even at this point, stochasticity in 
the earlier evolution can still influence the final equilibrium 
configuration. 

We would argue that more realistic initial conditions 
and physics (e.g. a mass spectrum of stars, a clumpy live gas 
potential, feedback, etc.) can only work to make the prob- 
lem more stochastic, rather than less. Therefore we contend 
that it is impossible to give a 'standard' outcome for any 
statistically identical clusters: some will be destroyed, while 
others will survive with different mass and structure. For 
this reason it is impossible to go backwards from an ob- 
served surviving cluster to its initial conditions except in a 
very broad statistical sense. Therefore we urge caution in de- 
riving the initial cluster population from an observed naked 
cluster population. 
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